{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This notebook produces Figure 1 in the main text, as well as Figure 1 in Further Supplemental Information."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Import Libraries and Set Project Directory"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Data Analysis\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import janitor\n",
    "import pandas_flavor as pf\n",
    "\n",
    "# Plotting\n",
    "import matplotlib.pyplot as plt \n",
    "import seaborn as sns; sns.set(color_codes=True)\n",
    "import matplotlib as mpl\n",
    "sns.set_style('ticks')\n",
    "sns.despine()\n",
    "\n",
    "# Solution to font issues:\n",
    "# https://stackoverflow.com/questions/52673444/\n",
    "# matplotlib-rcparams-not-recognizing-times-new-roman-mac-high-sierra\n",
    "plt.rcParams['font.family'] = 'serif'\n",
    "plt.rcParams['font.serif'] = ['Times New Roman'] + plt.rcParams['font.serif']\n",
    "\n",
    "# CHANGE THE FOLLOWING \"PROJECT_DIR\" FILEPATH TO REFLECT LOCATION OF REPLICATION FOLDER\n",
    "PROJECT_DIR = '/Users/jjares/Documents/research/Farm Subsidies/replication_materials'\n",
    "GRAPHICS_OUTDIR = f'{PROJECT_DIR}/graphics'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Produce Figure 1 (Main Text)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Read in data. For readability, subset to five main program groups."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "disaster_subcategories = [\n",
    "    'crop_disaster_program_of_2005_2007', 'elap', 'emergency_forestry_restoration_program',\n",
    "    'livestock_compensation_program_of_2005_2007', 'livestock_forage_disaster_program',\n",
    "    'livestock_indemnity_program', 'misc_disaster', 'non_insured_assistance_program',\n",
    "    'supplemental_revenue_assistance_program', 'tree_assistance_program', 'whip',\n",
    "    'emergency_conservation_program'\n",
    "]\n",
    "by_program_year = (\n",
    "    pd.read_csv(f'{PROJECT_DIR}/data/USDA Payment Totals by Program and Calendar Year (2004-2020).csv')\n",
    "      .assign(disaster = lambda d: d[disaster_subcategories].sum(axis=1),\n",
    "              arc_plc = lambda d: d['arc'] + d['plc'],\n",
    "              mfp = lambda d: d['2018_mfp'] + d['2019_mfp'])\n",
    "      .select_columns(['year', 'mfp', 'arc_plc', 'dcp_acre', \n",
    "                       'disaster', 'crp_annual_rental'])\n",
    "      .set_index('year')\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Save time series plot."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(10, 6))\n",
    "program_labels = {\n",
    "    'mfp': 'Market Facilitation Program',\n",
    "    'arc_plc': 'ARC/PLC',\n",
    "    'dcp_acre': 'DCP/ACRE',\n",
    "    'disaster': 'Disaster Programs',\n",
    "    'misc_commodity': 'Misc. Commodity Programs',\n",
    "    'crp_annual_rental': 'CRP Annual Rental',\n",
    "}\n",
    "program_markers = {\n",
    "    'mfp': 's',\n",
    "    'arc_plc': 'P',\n",
    "    'dcp_acre': 'o',\n",
    "    'disaster': '^',\n",
    "    'misc_commodity': '',\n",
    "    'crp_annual_rental': '*'    \n",
    "}\n",
    "\n",
    "# Plot vertical lines indicating implementation \n",
    "# of 2014 Farm Bill and start of trade conflict\n",
    "plt.axvline(x=2015, linestyle='--', color='black', alpha=0.5)\n",
    "plt.axvline(x=2018, linestyle='--', color='black', alpha=0.5)\n",
    "ax.annotate('2014 farm bill\\ntakes full effect', xy=(2014.95, 10), xytext=(2011.5, 9.9),fontsize=12,\n",
    "            arrowprops=dict(arrowstyle=\"->\", color='black'))\n",
    "ax.annotate('US-China trade\\nconflict begins', xy=(2017.9, 11.4), xytext=(2015.2, 12.6),fontsize=12,\n",
    "            arrowprops=dict(arrowstyle=\"->\", color='black'))    \n",
    "\n",
    "program_order_1 = ['dcp_acre', 'arc_plc', 'mfp', 'crp_annual_rental', 'disaster', 'misc_commodity']\n",
    "for i, program in enumerate(program_order_1):\n",
    "    if program == 'misc_commodity':\n",
    "        continue\n",
    "    plot = by_program_year.query(f'{program} > 10**7') / (10**9)\n",
    "    marker = program_markers[program]\n",
    "    if marker in ['*', '|']:\n",
    "        marker_size = 9\n",
    "    elif marker == 'P':\n",
    "        marker_size = 8\n",
    "    else:\n",
    "        marker_size = 6\n",
    "    sns.lineplot(x=plot.index, y=plot[program], \n",
    "                 label=program_labels[program], \n",
    "                 marker=marker,\n",
    "                 markeredgecolor=sns.color_palette()[i],\n",
    "                 markersize=marker_size,\n",
    "                 ax=ax)\n",
    "        \n",
    "ax.spines['right'].set_visible(False)\n",
    "ax.spines['top'].set_visible(False)\n",
    "plt.ylim(bottom=0)\n",
    "plt.ylabel('Billions of 2020 Dollars', fontsize=12)\n",
    "plt.xlabel('')\n",
    "plt.legend(loc='upper left', fontsize=12)\n",
    "plt.savefig(\n",
    "    f'{GRAPHICS_OUTDIR}/Figure 1.tiff', \n",
    "    dpi=1000, \n",
    "    bbox_inches='tight'\n",
    ")\n",
    "plt.savefig(\n",
    "    f'{GRAPHICS_OUTDIR}/Figure 1.pdf',\n",
    "    bbox_inches='tight'\n",
    ")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Produce Figure 1 (Further Supplemental Information)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "deflator = (pd.read_excel(f'{PROJECT_DIR}/data/CPI for All Urban Consumers (2000-2021, CUUR0000SA0).xls', \n",
    "                          skiprows=11)\n",
    "              .clean_names()\n",
    "              .assign(deflator = lambda d: d.set_index('year').loc[2020, 'annual'] / d['annual'])\n",
    "              .query('2004 <= year <= 2020')\n",
    "              .set_index('year')\n",
    "              .loc[:, 'deflator'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ers_program_recategorization = {\n",
    "    'Direct government payments': 'total_payments',  \n",
    "    'Production flexibility contract payments': 'direct',\n",
    "    'Fixed direct payments': 'direct', \n",
    "    'Cotton Transition Assistance Payments (CTAP)': 'direct',\n",
    "    'Cotton Ginning Cost-Share (CGCS) Program': 'direct', \n",
    "    'Average Crop Revenue Election Program (ACRE)': 'counter_cyclical',\n",
    "    'Price Loss Coverage (PLC)': 'counter_cyclical',\n",
    "    'Agriculture Risk Coverage (ARC)': 'counter_cyclical',\n",
    "    'Counter-cyclical payments': 'counter_cyclical',\n",
    "    'Loan deficiency payments': 'marketing_loan_alt',\n",
    "    'Marketing loan gains': 'marketing_loan_alt',\n",
    "    'Certificate exchange gains': 'marketing_loan_alt',\n",
    "    'Peanut quota buyout payments': 'direct',\n",
    "    'Milk income loss payments': 'ad_hoc_relief',\n",
    "    'Dairy Margin Coverage Program': 'ad_hoc_relief',\n",
    "    'Tobacco Transition Payment Program': 'direct',\n",
    "    'Conservation': 'conservation',\n",
    "    'Biomass Crop Assistance Program (BCAP)': 'direct',\n",
    "    'Supplemental and ad hoc disaster assistance': 'ad_hoc_relief',\n",
    "    'Market Facilitation Program': 'ad_hoc_relief',\n",
    "    'Miscellaneous programs': 'ad_hoc_relief'\n",
    "}\n",
    "ers_payments = (pd.read_excel(f'{PROJECT_DIR}/data/Gov Payments by Program (1933-2019, Nominal).xlsx', \n",
    "                             sheet_name='United States', skiprows=[0, 1, 2, 4], engine='openpyxl')\n",
    "                  .rename_column('United States', 'program')\n",
    "                  .remove_columns(['Unnamed: 4'])\n",
    "                  .set_index('program')\n",
    "                  .loc[list(ers_program_recategorization)]\n",
    "                  .fillna(0)\n",
    "                  .reset_index()\n",
    "                  .assign(category = lambda d: d['program'].map(ers_program_recategorization))\n",
    "                  .remove_columns(['program'])\n",
    "                  .groupby('category')\n",
    "                  .sum()\n",
    "                  .transpose()\n",
    "                  .rename_axis('year', axis=0)\n",
    "                  .rename_axis(None, axis=1)\n",
    "                  .reset_index()\n",
    "                  .assign(year = lambda d: d['year'].str.strip('F'))\n",
    "                  .change_type('year', 'int')\n",
    "                  .query('(2004 <= year) & (year <= 2020)')\n",
    "                  .assign(deflator = lambda d: d['year'].map(deflator)))\n",
    "\n",
    "# Verify that the program categories add up to the total in each year\n",
    "categories = [\n",
    "    'direct', 'counter_cyclical', 'ad_hoc_relief', \n",
    "    'conservation', 'marketing_loan_alt'\n",
    "]\n",
    "assert (ers_payments[categories].sum(axis=1) / ers_payments['total_payments']).between(0.999, 1.001).all()\n",
    "\n",
    "# Convert all program totals to 2020 dollars, and truncate at zero\n",
    "for field in categories + ['total_payments']:\n",
    "    ers_payments[field] = np.where(\n",
    "        ers_payments[field] >= 0,\n",
    "        ers_payments[field] * ers_payments['deflator'],\n",
    "        0\n",
    "    ) \n",
    "ers_payments = ers_payments.drop(columns='deflator')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "net_farm = (pd.read_excel(f'{PROJECT_DIR}/data/Net Farm Income (1910-2021F, Nominal).xlsx', \n",
    "                          skiprows=[0], engine='openpyxl', \n",
    "                          sheet_name='United States')\n",
    "              .rename_column('United States', 'field')\n",
    "              .remove_columns(['Unnamed: 4'])\n",
    "              .query('field == \"Net farm income\"')\n",
    "              .set_index('field')\n",
    "              .transpose()\n",
    "              .rename_axis('year')\n",
    "              .reset_index()\n",
    "              .assign(year = lambda d: d['year'].str.strip('F'))\n",
    "              .change_type('year', 'int')\n",
    "              .rename_column('Net farm income', 'net_farm')\n",
    "              .change_type('net_farm', 'int')\n",
    "              .query('(2004 <= year) & (year <= 2020)')\n",
    "              .assign(deflator = lambda d: d['year'].map(deflator),\n",
    "                      real_net_farm = lambda d: d['net_farm'] * d['deflator']))\n",
    "\n",
    "net_income = (net_farm.select_columns(['year', 'real_net_farm'])\n",
    "                      .merge(ers_payments, on='year')\n",
    "                      .eval('real_net_farm_before_payments = real_net_farm - total_payments')\n",
    "                      .eval('prop_net_farm = total_payments / real_net_farm'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_order = [\n",
    "    'real_net_farm_before_payments',\n",
    "    'direct', 'counter_cyclical', 'ad_hoc_relief', \n",
    "    'conservation', 'marketing_loan_alt'    \n",
    "]\n",
    "item_labels = [\n",
    "    'Net Farm Income Before Government Payments', \n",
    "    'Unconditional Direct Payments',\n",
    "    'Counter-Cyclical Commodity Payments',\n",
    "    'Ad Hoc Relief (e.g. Disaster, MFP, CARES)',\n",
    "    'Conservation',\n",
    "    'Marketing Loan Option Exercises'\n",
    "]\n",
    "plot_data = net_income.set_index('year')[plot_order].divide(10**6)\n",
    "fig, ax = plt.subplots(1, figsize=(7,7))\n",
    "#ax = net_income.set_index('year')[plot_order].divide(10**6).plot.area()\n",
    "ax.stackplot(\n",
    "    np.arange(2004, 2021), \n",
    "    [plot_data[field] for field in plot_order],\n",
    "    baseline='zero',\n",
    "    colors=[sns.color_palette()[i] for i in [0, 1, 4, 3, 2, 5]],\n",
    "    linewidth=0,\n",
    "    labels=item_labels\n",
    ")\n",
    "plt.ylabel('Billions of 2020 Dollars', fontsize=12)\n",
    "ax.spines['right'].set_visible(False)\n",
    "ax.spines['top'].set_visible(False)\n",
    "ax.collections[0].set_alpha(0.5)\n",
    "ax.legend(loc=(0.1, 0.05))\n",
    "plt.savefig(f'{GRAPHICS_OUTDIR}/Further Supplemental Information Figure 1.pdf', bbox_inches='tight')\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
